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Tomography— cross-sectional imaging based on measuring radiation transmitted through an 
object along different directions— enables non-invasive imaging of hidden stationary objects, 
such as internal bodily organs, from their sequentially measured projections. Here we adapt 
tomographic methods to visualize— in one laser shot— the instantaneous structure and 
evolution of a laser-induced object propagating through a transparent Kerr medium. We 
reconstruct 'movies' of a laser pulse's diffraction, self-focusing and filamentation from phase 
'streaks' imprinted onto probe pulses that cross the main pulse's path simultaneously at 
different angles. Multiple probes are generated and detected compactly and simply, making 
the system robust, easy to align and adaptable to many problems. Our technique could 
potentially visualize, for example, plasma wakefield accelerators, optical rogue waves or fast 
ignitor pulses, light-velocity objects, whose detailed space-time dynamics are known only 
through intensive computer simulations. 
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n the 1870s, English photographer Eadweard Muybridge used 
an array of cameras activated by trip lines strung across a 
I horse's path to show that it became temporarily airborne while 
galloping, thus settling a hotly debated question of his time 1 . 
Muybridge's demonstration did not rely on repeatability of the 
horse's motion. Rather, he established his point by capturing a 
series of stop-action frames within one cycle of the horse's gallop, 
a few of which revealed the horse's legs collected beneath its 
airborne body. In modern science, researchers still seek a method 
to record a Muybridge -like multiframe movie of evolving light- 
velocity objects created by laser or charged particle pulses during 
a single transit through matter. These objects include optical 
rogue waves used to model giant ocean waves 2 , nonlinear plasma 
wakes used to accelerate charged particles 3 , fast ignitor pulses 
used to trigger laser fusion 4 , filamenting laser pulses used for 
atmospheric analysis 5 and filamenting electron bunches used to 
model astrophysical jets 6 . Sophisticated single-shot techniques 
exist for measuring the complete electric field profile 
E (z(loc) )X (±)) (or charge-density profile p(z (loc) ,x (1) )) of an 
ultrashort light 7-9 (electron 10 ) pulse that initiates such events, 
analogous to the opening frame of Muybridge's movie. Here, 
^(loc) _ z _ V f ^ s foe local longitudinal position relative to the 
centre of the pulse propagating at velocity v, a vector 

transverse to the propagation (z) axis, z distance into the 
medium and t time after the pulse enters the medium. However, 
when the pulse propagates nonlinearly through matter, its profile 
E(z (loc) ,x (1) ,t) (or p(z (loc) ,x (1) ,0)— and the co-propagating 
refractive index profile An(z^ oc \x^\t) that it creates in the 
medium by, for example, Kerr effect, ionization, or plasma wave 
generation — evolve in complex, often unrepeatable, ways for 
which no single-shot recording method exists. Detailed 
knowledge of such evolving structures is essential for 
understanding, optimizing and scaling the myriad applications 
of laser- and particle-matter interactions, yet is available for 
individual shots only from computer simulations based on 
estimated initial conditions. 

Several methods to visualize laser- and particle-generated 
structures were developed previously. For example, 'snapshot' 
images of quasi-static (^-independent) laser-generated structures 
An(z^ oc \x^^) (refs 11,12) — analogous to Muybridge's horse 
standing on a light-velocity treadmill— were produced using 
frequency- domain holography (FDH) 13,14 , in which a 
longitudinally stretched probe pulse (and a reference pulse with 
which it interfered in a spectrometer) copropagated with the 
object. Such images revealed relativistic curvature of quasi-linear 
plasma wake fronts 11 and resolved molecular and plasma 
contributions to formation of a stable laser filament in gas 12 . 
However, FDH snapshots blur if the objects evolve as they 
propagate. Thus, strongly nonlinear plasma wakes 15 yielded FDH 
snapshots 16 that averaged over time variations essential to their 
function 17 ; similarly, FDH cannot image the stochastic break-up 
of a pulse into multiple filaments. A frequency- domain streak 
camera (FDSC), in which the evolving object imprints a phase 
streak on an FDH-like probe pulse that crossed its path 
obliquely 18 , yielded partial single-shot information about the 
object's evolution, but not a Muybridge-like image sequence. 
Motion-picture images of time- varying light- velocity objects have 
therefore been produced only by probing them transversely at 
different time delays over multiple shots 19-21 , analogous to 
strobing Muybridge's horse over successive gallop cycles. One 
limit of this approach for light- velocity objects is that the transit 
time t of a transverse probe pulse across the object's width 
averages internal longitudinal structure of scale Az^ oc ^<vt. 
Nevertheless, this approach has achieved impressive z^ oc ^ 
resolution of very thin objects by using extremely short probe 
pulses. For example, Buck et al. 20 produced a shadowgraphic 



movie of plasma wakes of period Az( loc ) < 4 urn and width 
~ 5 Jim by probing transversely with ~ 6 fs pulses on successive 
shots. Nevertheless, the principal limit of multishot visualization 
is that it is not viable for non-repetitive or stochastic events 
typical of highly nonlinear interactions 2 ' 15 ' 21 , propagation in 
turbulent media , or interactions driven by sources with low 
repetition rate 23 , unstable pointing, or other shot-to-shot 
fluctuations. 

Here we present a new approach to visualizing evolving light- 
velocity objects that yields a Muybridge-like movie sequence in 
one shot, enabling up to four- dimensional (4D) visualization of 
infrequent or stochastic laser- or particle-driven events for the 
first time, and greatly speeding visualization of events that are 
stable from shot to shot. The method utilizes multiple probe 
pulses, each an FDSC 18 crossing the object's path at a different 
angle. However, the heart of the method is based on two key 
innovations that go well beyond FDSC. First, to avoid runaway 
complexity and cost, we developed compact, inexpensive methods 
to generate and detect multiple synchronized probe pulses. 
Instead of an unwieldy hyper-Michelson interferometer , we 
generated the multiprobe array in a single step by cascaded four- 
wave mixing. On the detection side, instead of multiple expensive 
spectrometers, we multiplexed all probes to one spectrometer and 
extracted all phase streaks in one step from one high- 
information -density spectral hologram. Consequently, our 
apparatus is as simple and robust as a standard pump-probe 
experiment. Second, to transform phase streaks into a movie, we 
adapted established methods of computerized tomography 25 ' 26 
for the first time to reconstruct ultrafast spatiotemporal 
dynamics. We generalized the Radon transformation — 
conventionally used to reconstruct the spatial structure of a 
stationary object (for example, an internal bodily organ, 
underground petroleum reserve, or hidden manufacturing 
defect) from its projections, measured sequentially 25 or 
simultaneously 27 — to reconstruct both the spatial structure and 
the ultrafast temporal dynamics of an evolving light-velocity 
object. Similar to any optical visualization system, frequency- 
domain tomography (FDT) has resolution limits in z^ oc \ x^ ^ 
and t. Instead of transit time, z^ oc ^ resolution is limited by probe 
bandwidth, which can be increased as needed by supercontinuum 
generation 14 . As for a simple camera, resolution is limited 
by diffraction and lens aberrations. Finally, evolution time 
t resolution is limited to changes occurring over propagation 
distances z larger than the object's dimensions, analogous to 
the paraxial approximation embedded in most mathematical 
treatments of nonlinear pulse propagation, yet allowing 
visualization of most dynamics of interest. 

Results 

Phase streaks and projection angles. Figure la schematically 
depicts a single probe pulse of centre wavelength 2 pr crossing the 
path of an evolving pump -generated object An(z^ c \ x oh , y ohy t oh ) 
at laboratory angle 0, thereby accumulating phase streak 

> x pr> >pr) via cross-phase modulation. Here An and \jj 
are expressed in the coordinates of the respective co -moving 
frames of object (ob) and probe (pr) pulses. We shall assume that 
the effective propagation length L eff over which any region of the 
probe profile overlaps the object is shorter than a probe diffrac- 
tion length Lfifi = n(Axl h + Ayl h )/X pr n(X pr ), which is valid for the 
experiments presented below, and for a wide range of other 
interactions of interest (for example, laser-plasma accelerators in 
gas jets). Here, Ax ob , Ay oh are the object's transverse radii, and L eff 
is less than the medium length L because the object drifts across 
the probe as they propagate at different angles and/or velocities. 
Diffraction of the probe pulse from the object can then be 
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Figure 1 | Phase streak and projection angle, (a) Schematic of chirped 
probe pulse crossing path of evolving index object An at laboratory angle 6, 
producing phase streak at projection angle </>. At f o b = 0, the index object 
(red circle) overlaps with the right side of the probe profile and the probe 
keeps going until t ob = L ob /v ob when it overlaps with the left side of the 
probe, leaving a phase streak across the probe, (b) Plot of </> versus 6 for 
three different probe velocities, as indicated in legend. The object 
propagates at v ob = 0.68c in all three cases. 

neglected, and its phase shift can be expressed as an integral 

lk {8) , W P r)= T / An(0) ( Z ob C) . *ob, Job, fob)Vp r df pr 

^pr JO 

(i) 

over the index object. In this limit, a point of the probe pulse 
accumulates a phase shift only when a linear chord of the object 
sweeps through it. Equation (1) is a 4D generalization of the 
Radon transform widely used in conventional tomography 25 . If 
diffraction is important, then equation (1) should be replaced 
with a Fresnel diffraction integral. 

To use equation (1), An(z^ ,x oh ,y oh , t oh ) must be transformed 
from the coordinates of its own co- moving frame to those of the 
probe frame (see Methods). The Radon transform equation (1) 
then becomes: 

^)(4^), Wpr )=^ / An z^cosfl + vsine 

^pr Jo L 

- Oob - v pr cos#)t ob , - z^smO 
+ x pr cos# - v w t oh sm6,y pY , t oh ]dt oh . (2) 



The phase streak (equation (2)) represents a projection of the 
index object in 4D space (z^ c \ x oh , y oh , t oh ) onto the three- 
dimensional (3D) probe profile space C4r° C ^*pr>7pr)- 

In conventional tomography, the projection angle is simply the 
laboratory angle 9 at which probe light intersects the stationary 
object. In FDT, the projection angle (j) (9) differs from, but 
depends on, 9. (j) is the angle between the object's propagation 
direction and the axis of the phase streak (see Fig. la). 
A straightforward derivation (see Methods) yields 

0= arg[(v ob - v pr cos#) + f(v pr v o b/c)sin0] , (3) 

where the argument in brackets is a complex number of the form 
x+iy corresponding to angle 4> — tan - (y/x) with the real axis. 
Figure lb plots 9 versus (j) at different v pr for fixed v ob . As 
typically v 0 b~v pr for light-velocity index objects, a relatively 
narrow range of 9, which is convenient for laboratory set-up, can 
yield a broad range of 0, which optimizes image reconstruction 
quality. For example, for v Q b = 0.68c and v pr = 0.66c, typical of an 
800-nm pump and 400-nm probe in glass, equation (3) yields 
projection angle range — 70° < (j)< 70° for probe angle range 
- 10°<fl<10°. 

Projection angles 0 = 0 and n/2 have special meaning. A phase 
streak with 0 = 0 is obtained with 9 — 0 (co-propagating object 
and probe) and v ob > v pr (equivalently, (ft = it corresponds to 
0 = 0, v 0 b<v pr ). Such streaks reveal evolution of the object's 
transverse (x^ 1 ^) profile, as it drifts longitudinally along the probe 
profile. Unlike FDH where velocity walk- off between co- 
propagating object and probe is a disadvantage because it blurs 
the image, in FDT it becomes an advantage because it creates a 
phase streak at 0 — 0 that is equivalent to streaks created by 
angular walk- off, and adds to the information pool for 
tomographic reconstruction. Streaks at 0~7i/2 reveal evolution 
of the object's longitudinal (z^ loc ^) profile as it drifts sideways 
across the probe profile. Optimal FDT reconstruction should 
include projection angles near both of these values. 

FDT movie reconstruction. Each probe pulse carries a 3D 
projection (equation (2)) of the complete 4D object that is, in 
principle, fully recoverable. Below, however, we present recording 
procedures that select a two-dimensional (2D) slice of the phase 
streak in equation (2) at a desired constant y pr —yo- Thus, from 
here on, we outline the reduced procedure for recovering a 
tomographic movie of the 3D slice An(z^ c \ x 0 \„ y 0 , £ Q b) from the 
selected 2D projections \jj^\z^r C \x w ,y 0 ) rather than the full 4D 
visualization. 4D imaging can be achieved by splitting the probes 
and sending copies to different spectrometers, each slicing a 
different y 0 . 

Of various conventional tomographic algorithms 25 , here we 
adapt algebraic reconstruction techniques (ARTs) 25,28 ' 29 , because 
they are compatible with a smaller number of probes, are 
less susceptible to reconstruction artefacts and more readily 
incorporate independent knowledge of the object than other 
Radon-transform-based algorithms. ART is an iterative algorithm 
in which one guesses the object's shape at each step, compares 
with measured projections, then modifies the guess until residuals 
are minimized (see Methods). In each iteration, independent 
information — for example, measured laser pulse profiles entering 
and exiting the medium and known optical properties of the 
medium — were incorporated and constrained trial solutions for 
An. No prior assumptions about the object's symmetry are made, 
in contrast to, for example, Abel inversions that assume 
cylindrical symmetry. 

FDT resolution limits. The ability of FDT to resolve an evolving 
object An(z^°\ x 0 \>, y 0 , £ Q b) must be considered separately for each 
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variable z^ , x Q b and t ob . Resolution of 'intraframe' (z^ } , x ob ) 
structure at a given time t oh is determined by imaging optics 
and light-wave properties, as for conventional imaging of 
non-evolving objects. Wide probe bandwidth optimizes z^ 
resolution. For example, to resolve non-evolving ~ 2 urn 
^ob^ s t ructures by FDH, supercontinuum probes were used 14 . 
Diffraction and lens aberrations determine x oh resolution of each 
frame. Single-element imaging lenses achieve ~ 5 jim resolution 
routinely. Aberration -corrected optics enable higher x oh 
resolution. 

The new criterion that enters with tomographic movies is 
'interframe' time (t oh ) resolution. FDT is designed to fully resolve 
objects that satisfy the paraxial approximation, that is, the drive 
pulse and its index object evolve slowly over propagation 
distances z 0 ^ of the order of the object's dimensions. Thus, for 
an object of tens of microns longitudinal and lateral extent 
(typical of those created by a loosely focused ~ 100 fs laser pulse), 
picosecond- scale evolution is fully resolved. More rapid evolution 
is still observable, with limited resolution. Interframe resolution 
depends on the length of the phase streak relative to the object's 
dimensions. FDT can resolve ~ N stages of the object's temporal 
evolution, where N is the number of separated objects of 
dimensions (z^ c \ Ax G b) that can be lined up along the streak 
axis. Technical details and estimates of interframe resolution for 
experiments discussed below are presented in Methods. 

FDT image quality and reconstruction artefacts. Image quality 
depends on the number and distribution of probes, and is 
quantified via the normalized root mean squared (RMS) deviation 
d (or 'error') of a reconstruction An(z^ oc \x, t) from the evolving 
object An 0 (z^ oc \x, t). We define d by generalizing its conven- 
tional 2D expression 25 ' 26 to 3D: 



Jloc) 



EgOocW [Mz (l0C U, t) - An 0 (z^\x, t)} 2 
\ E,aoc), x ,JAno(^),x,0-An 0 ] 2 



(4) 



The z^ oc \x sums run over intraframe pixels, the t sums over 
frames and An 0 is the object's average index. The actual object 
An 0 (z^ oc \ x, t) is, of course, not directly accessible in a real 
tomographic reconstruction. Direct comparison of the re- 
constructed object with physical simulations and independent 
measurements are thus important for validating FDT 
reconstructions. In addition, d can be determined for a given 
probe configuration by simulating the reconstruction of a 
'phantom' object that is known exactly. 

Phantom simulations. Figure 2a shows selected 2D (x oh versus 
z^) snapshots of an artificial phantom index object at seven 
different times t ob (listed along the top) after entering a medium. 
These times correspond to object centre positions ranging from 
entrance (z 0 \> = 0) to exit (z Q b = 3 mm) of the medium. The 
horizontal spatial scale of each snapshot denotes z^ c \ with the 
object's leading edge (analogous to the head of Muybridge's 
horse) to the right along the propagation direction z ob . The 
phantom does not evolve by a real physical process, although 
some of its general features (for example, propagation length, x ob 
and z^ dimensions, evolution speed) were chosen to resemble 
those that occur in experiments below. The phantom's detailed 
features were chosen to illustrate resolution limits, and to evaluate 
reconstruction artefacts, more effectively than a real physical 
process. Specifically, the object starts as a hollow rectangle with 
thin boundaries of widths Azq 0 ^ = 2 jim (left and right) and 
Ax 0 = 5 jim (top and bottom), as shown by dotted curves in 
Fig. 3b,c, respectively. Such a thin rectangle separately and 
stringently tests transverse and longitudinal resolution limits. 



As it propagates at v Q b = 0.68c, chosen to equal pump group 
velocity in experiments, the rectangle narrows along x ob over 
0<£ o b<7-3ps, thus mimicking self-focusing observed in those 
experiments. During a short transition period (£ ob ~7.3 ± 0.2 ps), 
a 'dot' of Gaussian profile appears to the lower left of the rec- 
tangle, quickly grows to An max and falls to 0.2An max (Fig. 3a, 
dotted curve), mimicking the time scale of plasma generation and 
partial recombination in the experiments. As the dot grows 
within an interval (~ 100 fs) comparable to the object's duration 
Az^ c Vv ob , it tests interframe resolution. By breaking axial sym- 
metry, it also tests the algorithm's ability to reconstruct objects 
without prior assumptions about symmetry. The narrowed 
rectangle then expands longitudinally over 7.3 ps< t ob < Hps, 
mimicking group -velocity stretching of a laser pulse. 

Figure 2b-d shows tomographic reconstructions of the 
phantom from three probe configurations I through III. 
Figure 3a-c compares selected lineouts of the reconstructed 
images, with those of the object (dotted curves), whereas Fig. 3d 
plots RMS errors versus iteration number. In all cases, the 
probe bandwidth and imaging optics fully resolved the 
object's intraframe (z^ oc \x) features. Variations among recon- 
structions I— III thus arose solely from different probe number, 
angular distribution and velocity as follows: (I) 19 probes 
each with v pr = 0.70c distributed equally over projection angles 
- 90° < 0 < 90°; (II) 18 probes each with v pr = 0.681c distributed 
equally over - 90° < 0 < 80°; (III) 5 probes each with v pr = 0.66c 
distributed equally over — 70° < (j) < 70°. As expected, configura- 
tions I, II with 18-19 probes yielded sharper frame images of the 
slowly evolving rectangle (Fig. 2b,c) and sharper lineouts of its 
edges (blue and red curves, respectively, in Fig. 3b,c), than the 
5 -probe configuration III, which yielded more blurred images 
with reconstruction artefacts (Fig. 2d) and less distinct edge 
lineouts (dashed green curves, Fig. 3b,c). Reconstructions I and II 
also converged more rapidly towards smaller steady- state RMS 
error than configuration III (Fig. 3d). Thus, large probe number is 
one factor that promotes high-fidelity reconstruction of slowly 
evolving objects. On the other hand, all three configurations 
resolved the rectangle's slow narrowing and lengthening equally 
well (see Fig. 2b-d). 

Reconstructions of the 'dot' reveal more subtle comparisons. 
Configuration I, having the widest (j) range and largest 
| v pr — v Q b| = 0.02c, best resolved the dot's ultrafast evolution 
(Fig. 3a, blue curve), consistent with the discussion of FDT 
resolution above and in Methods. Configuration III, having the 
same |v pr — v ob | and similar (j) range, resolved this feature nearly 
as well (Fig. 3a, dashed green curve), despite only five probes. 
Configuration II, in contrast, failed completely to resolve it 
(Fig. 3a, red curve), a consequence of its small |v pr — v ob | = 0.001c. 
Nevertheless, configuration II yielded sharper images of the 
rectangle (Fig. 2c) and smaller RMS error (Fig. 3d) than 
configuration I, despite one less probe. Evidently reduced blurring 
of slowly evolving objects is a compensating advantage of small 
I Vpr — v Q b|. Thus, in choosing |v pr — v ob |, an FDT system designer 
must compromise between resolving interframe evolution and 
minimizing intraframe RMS error. 



Experimental FDT setup and procedure. As a laboratory 
demonstration of FDT, we visualized a 3D slice An(z^ c \ x oby y 0 , 
t oh ) of the evolving nonlinear refractive index envelope of a pump 
pulse (duration T pu =100fs, wavelength ^ pu = 800nm, energy 
0.4<£ pu <0.7 uj) focused with /-number ~30 to radius 
w 0 — 25 jim (peak incident intensity 0.4 < I 0 < 0.7 TW cm ~ 2 ) 
near the entrance of a Kerr medium: a fused silica plate of 
linear index n 0 (l^ u ) = 1.45 and lowest-order Kerr coefficient 
n 2 ~2x 10 ~ 16 cm W -1 . Incident power exceeded the critical 
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Figure 2 | Phantom simulations, (a) 2D snapshots (x ob versus z^) of evolving 'phantom' object at seven selected times t ob after entering a medium 
(that is, seven locations z ob = v ob f ob along its propagation path) as it propagates with velocity v ob = 0.68c. The longitudinal axis of each snapshot 
indicates local distance from the centre of the object, with the object's leading edge to the right. Remaining rows show tomographic reconstructions 
of the phantom object for three probe configurations: (b) Nineteen probes at projection angles -90°<(/x90° with 10° separation, v pr = 0.70c 
(configuration I); (c) 18 probes at -90°<(/x80° with 10° separation, v pr = °- 681c (configuration II); (d) 5 probes at -70°<(/x70° with 35° 
separation, v pr = 0.66c (configuration III). The colour bar shows the dimensionless refractive index change of original and reconstructed objects. 
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Figure 3 | Analysis of phantom simulations, (a) Peak An of the Gaussian 'dot' versus evolution time £ 0 b = z ob A ob , demonstrating 
interframe resolution; (b,c) lineouts at f ob = 7.3ps (that is, z ob = 1.5mm) of the rectangle along z o '^ (b) and x ob (c), demonstrating intraframe 
resolution. Curves in a-c refer to phantom simulations in Fig. 2: original object (black dotted curve) and reconstructions with 19 (blue solid curve), 
18 (red dashed curve) and 5 (green dash-dot curve) probes, (d) Normalized root mean square error of the same three tomographic reconstructions 
versus iteration number. 
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power P cr = 3.7/lp U /87in 0 n 2 ~3.2 MW for self-focusing by a factor 
of 1.2 to 2.2. Our plate thickness L = 3mm satisfied the weak 
diffraction criterion L e ff<Lcijf? underlying equation (1). At our 
highest £ pu , we observe a 3D slice of dynamics preceding self- 
guided filament formation 30 ' 31 . 

The primary technical challenge in implementing FDT in the 
laboratory is to avoid runaway complexity and cost in generating, 
formatting and detecting a multiprobe pulse array. Conventional 
multiprobe experiments require an array of beam splitters to 
divide probes from pump, and a multimirror hyper-Michelson 
interferometer to format the probe train. Such setups are 
challenging to align and sensitive to vibrations. Moreover, 
recording the phase streaks by conventional FDH methods 
requires a separate spectrometer with charge-coupled device 
(CCD) detector for each probe. 

We addressed this challenge with the setup in Fig. 4. Only two 
probe-generating' pulses (800 nm, 30 fs, 30 uj) were split directly 
from the pump. These crossed simultaneously at a small 
adjustable angle (a ~ 5 mrad) in a three-layer structure consisting 
of a P-barium borate (BBO) crystal sandwiched between two 
HZF4 glass plates. Cascaded four-wave mixing in the first HZF4 
plate (5 mm thick, 4 HZF4) ~ io- 15 cm 2 W l ) created a fan of up 
to eight 800 nm daughter pulses. The number depended on probe 
generator intensity and was limited by the onset of self-focusing 
and self-phase modulation. The BBO crystal (Type I, 500 Jim 
thick) then frequency- doubled each beam and created an 
additional probe midway between each fundamental pair by 
sum-frequency generation, producing a fan of up to fifteen 
400 nm probes separated from each other by angle a/2. The 
second HZF4 plate (15 mm thick) chirped them to 600 fs 
duration. Any desired subset of the array was easily selected by 
blocking unwanted probes. Here, for reasons discussed in 
Methods, we selected five probes at 0 = 0.1°, 1.4°, -1.2°, 



— 7.6° and 9.5° inside the fused silica, corresponding to 
projection angles 0 = 1.0°, 27°, -25°, -65° and 68°, respec- 
tively, for v G b = 0.68c and v pr = 0.66c. 

Spatiotemporal overlap of multiple probes was achieved 
automatically by imaging the first HZF4 glass plate to the sample 
with lens L\ (f=20, 3 cm diameter). The pump -generated object 
swept across and imprinted a phase streak on each probe. To 
probe at angles 9 beyond the aperture of L u mirrors can re- direct 
selected probes from the fan along independent delivery lines to 
the sample. This was done for our two largest angle probes P 4 , P 5 . 
As shown in Fig. 4, a chirped 400-nm reference pulse also split 
from the pump, co-propagated with it, advanced temporally by 
T~2-4ps. 

Lens L 2 (f— 15 cm, /-number 5.6) imaged reference and phase- 
modulated probes from the sample exit face to the slit of a single- 
imaging spectrometer, which selected a lineout of constant y 0 , 
thus limiting the final reconstructed objects to 3D slices 
An(z^ c \x oh ,y 0) t oh ) of the full 4D object An(z^ c \x ohy y oh , t oh ). 
Lens L 2 also ensured simultaneous, spatially overlapped delivery 
of all probes within its aperture. Reference and probes interfered 
inside the spectrometer, projecting a grid-like frequency- domain 
intensity pattern or hologram I(cd,x) (Fig. 4b) onto the CCD. This 
hologram stored the phase modulations of all probes in one shot, 
analogous to coherent multiplexing methods in holographic data 
storage 32 . To analyse phase shift in each probe, a 2D Fourier 
transform of I(co,x) yielded a reciprocal 2D hologram J(T,0,.) 
(Fig. 4c), in which each reference-probe interference pattern 
appeared as an isolated peak (P 1 through P 5 ) at a position 
determined by probe angle 6 i and time- delay T. Peaks P 4>5 appear 
at a different T than other peaks because P 4 and P 5 have their own 
references (not shown). Each phase streak was reconstructed by 
windowing and inverse Fourier-transforming a specific peak P b as 
in FDH. The peak P { must be well separated to avoid cross-talk. 



Spectrometer 




Temporal delay (ps) 

Figure 4 | Schematic set-up for single-shot FDT of evolving laser-generated structures. An 800-nm pump pulse of intensity /, coupled into and out of a 
glass slab by dichroic mirrors, creates an evolving luminal-velocity refractive index structure due to the nonlinear refractive index n 2 l and pump-generated 
plasma. Two probe-generating pulses split from the pump generate an array of chirped, frequency-doubled probe pulses in an HZF4/BBO/HZF4 
sandwich. Lens L-| directs each probe across the pump path through the glass, where each acquires a phase streak that records evolution of the refractive index 
structure from a different angle. Optionally, selected probes can be diverted to angles beyond the optical aperture of (see probes P 4; P 5 ), or blocked as 
shown (see dashed optical paths). Lens L 2 images undeflected, unblocked probes PrP 3 from the glass exit face to the entrance slit of an imaging 
spectrometer, which selects a slice of the projected image. Large-angle probes P 4 , P 5 were relayed by confocal lenses equivalent to and L 2 (not shown). The 
probes interfere inside the spectrometer with a temporally advanced 400 nm reference pulse, creating a 2D frequency domain hologram on a CCD at the 
spectrometer's detection plane, shown in b for five probes. Fourier transformation of this data yields a reciprocal hologram, c, in which peaks P] through 
P 5 encode the phase modulations of the five probes. Phase streaks in Fig. 5 are recovered by windowing and inverse Fourier transforming these five peaks. 
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This requirement, together with the angular acceptance of the 
CCD and imaging optics, governed our selection of five probes 
(see Methods). Although this five-probe system adequately 
resolved the main self-focusing dynamics, a substantial increase 
in probe number and image quality is available by investing in 
large aperture, aberration- corrected lenses and a larger, more 
finely pixelated CCD. 

Figure 5 shows five phase streaks for £ pu = 0.7uJ. Small (j) 
streaks (top row) highlight evolution of the object's transverse 
(*ob) profile. Oscillations in transverse radius Ax oh and peak 
index change An max , indicating dynamic balance between self- 
focusing and defocusing, were evident in the last ~ 1/3 of these 
streaks. Streaks near 6 ~ 70° (bottom row) highlight evolution of 
its longitudinal (z^ ) profile, which remained nearly constant 
due to the low dispersion of fused silica. To visualize evolution of 
the object's full An{z^ c \x oh ,y 0 ,t oh ) profile, we tomographically 
reconstructed a movie from all five streaks using the ART 
algorithm 25 ' 28 ' 29 . 



Experimental movies of nonlinear laser propagation in glass. 

Figure 6a shows movie frames, or 2D snapshots, of the nonlinear 
index profile An(z^ c \ x G b) at five selected propagation times £ Q b 
after entering (£ ob =0) and before exiting (£ ob = 15 ps) the Kerr 
medium, for E pu from 0.4 uj (top row) to 0.7 uj (bottom row). 
The interframe spacing (At ob = 2.4ps, Az ob = 500 jim) approxi- 
mates the interframe resolution limit for transverse profile 
variations and twice this limit for longitudinal variations. 
Supplementary Movie 1 shows continuously streaming movies 
with higher frame density. The main feature in each frame 
is a positive An(z^ c \x 0 b) profile of transverse lie 2 radius 
13 < Ax ob <25 Jim and longitudinal duration ~ 20 um, or 
~ 100 fs, that is attributable mainly to the instantaneous lowest- 
order nonlinear Kerr response n 2 I pu of fused silica to the pump 
pulse. Diffraction, characterized by length L^g =7i(Ax ob ) / 
^p U n 0 (^ pu ), and self- focusing, characterized by focal length 
L n i — ^pJ2nn 0 n 2 Ip U (ref. 33), respectively, govern most transverse 
pump dynamics. The reconstruction resolves them as long as 
L^>500|im (that is, Ax ob >14um) and L nl >500jim (that 



is, I pu < 1 TWcm -2 ). Self-focusing beyond these limits can 
introduce dynamics faster than the interframe resolution, as 
well as additional nonlinearities such as plasma generation 12 ' 34 ' 35 
and higher-order Kerr effect 36-38 . Dispersion, characterized 
by length L dis = t^J fi 2 — 277 mm, governs evolution of the longi- 
tudinal profile 33 . As L dis »L, the pulse and its n 2 I pu profile 
propagate with negligible change in duration. 

The reconstructed movies in Fig. 6a depict various propagation 
regimes. For £ pu = 0.4uJ (top row), Ax oh and index shift 
An(Zjf^ — 0, x ob — 0) at the centre of the profile remained nearly 
constant over 15 ps (3 mm propagation distance), indicating that 
diffraction and self-focusing were balanced. For £ pu = 0.5uJ 
(second row), Ax ob contracted slightly, whereas An(0, 0) increased 
slightly (J pu increased from 0.5 to 0.6 TWcm -2 ), indicating that 
self-focusing slightly dominated. In these cases, the exit pump 
spectrum (Fig. 6b, top two rows) and spatial profile (Fig. 6c, top 
two rows) retained nearly their incident shapes (not shown), and 
transverse dynamics were fully resolved. At E pu = 0.6 uj, the pulse 
self-focused strongly and monotonically to Ax ob ^15jim and 
I pu > 1.5 TWcm , indicating that self-focusing dynamics 
slightly exceeded the interframe resolution limit near the end of 
the medium. New structure developed in the output spectrum 
(Fig. 6b, third row), a consequence of self-phase modulation, and 
in the spatial profile, which split into three lobes (Fig. 6c, third 
row). The last movie frame shows the uppermost of these lobes 
separating from the central profile at t oh = 12.2 ps (z ob = 2.5 mm), 
thus capturing the onset of multifilamentation. The other lobe lies 
outside the image plane. 

At £ pu = 0.7uJ, similar structure developed in the exit 
spectrum and beam profile (Fig. 6b,c, fourth row). In this 
case, however, the pulse self-focused to Ax ob ^15jim 
(J pu > 1.5 TWcm" 2 ) within 7.4 ps (z ob =1.5 mm; Fig. 6a, fourth 
row) instead of 15 ps (z ob = 3mm), after which further collapse 
was arrested and substructure developed in the index profiles. 
One such substructure is the split- off of the upper spatial lobe, 
now evident at f ob = 9.8ps (z ob = 2.0mm). The dominant new 
feature, appearing at £ ob =12.2ps (z ob = 2.5mm), is a steep- 
walled index 'hole' near the centre of the An profile. No such 
'hole' was observed in the exit beam profile (Fig. 6c, fourth row). 



0 = 0.1° 0=1.4° 0 = -1.2° 




-80 -40 0 40 80 -80 -40 0 40 80 -80 -40 0 40 80 




-60 -30 0 30 60 -60 -30 0 30 60 



(loc) / \ (loc) / x 

Figure 5 | Phase streaks induced by the evolving index profile. E pu = 0.7uJ and projection angles are 0 = 1.0°, 27°, -25° (top row), and -65° and 
68° (bottom). Vertical (horizontal) scales denote transverse (longitudinal) position of the pump-induced phase streak within the temporally stretched 
probe pulse profile. The right-most end of each streak corresponds to the entrance of the medium; x pr = Zp r 0C ^ =0 is approximately the midpoint of the 
pump pulse propagation through the Kerr medium. The spectrometer slit was centred on the images that L 2 projected at the spectrometer entrance. 
Colour bars give phase shift in rad. 
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Figure 6 | Single-shot tomographic movies of the evolving index profile, (a) Selected 2D snapshots An(z { ob \x ob ) of the pump nonlinear index profile 
at five different propagation times f ob indicated at the top and pump energies indicated at left. The horizontal axis of each snapshot is given in local 
position z^ c \ The colour bar shows the dimensionless refractive index change, (b) Spectra of the transmitted pump pulses. The spectrum in the top row is 
nearly identical to the incident spectrum, (c) Near-field images of the transmitted pump spatial profiles, with a resolution of ~20|im. Scale bar, 30 |im. 



Thus, it is the result of a negative index change that locally cancels 
the positive nonlinear index change n 2 I pu . This feature appears 
within a single frame and thus, similar to the 'dot' in the phantom 
simulations of Fig. 2, tests the interframe resolution of the 
tomographic reconstruction. 



Discussion 

To help understand the results in Fig. 6, we calculated 
propagation of the pump using a nonlinear Schrodinger 
equation (NLSE) 30 ' 31 that included diffraction, group -velocity 
dispersion, self- steepening, electronic and Raman-induced Kerr 
self-focusing, multi-photon and tunnelling ionization 39 , electron- 
hole recombination and plasma defocusing (see Methods). Input 
pulses were modelled as Gaussians that retained cylindrically 
symmetry as they propagated 30 ' 31 . The calculations thus do not 
capture asymmetric multifilamentation observed in some FDT 
images, but adequately model propagation up to the formation of 
such filaments. 

Figure 7a shows multiframe movies of the evolving An(z^ c \ 
x ob> t 0 b) from NLSE simulations neglecting high-order Kerr effect, 
for direct comparison with FDT reconstructions in Fig. 6a. 
Simulated and reconstructed movies agree with high fidelity 
despite the lack of any adjustable parameters. To quantify this 
agreement, Table 1 lists RMS deviations d between the FDT 
reconstructions in Fig. 6a and the NLSE simulations in Fig. 7a. 
For most frames of Fig. 6a, d<0.5, even smaller than for the 18- 
19 probe configurations of the phantom simulation (see Fig. 3d). 
This reflects the high fidelity of these reconstructions, as well as 
the absence of sharp -edged structures that were inserted into the 
phantom simulations to test stringently for reconstruction 
artefacts. We observe d>0.5 in two places. First, frames that 
encompass the entrance (£ ob <2.5ps, z ob <0.5cm) or exit 
(£ 0 b > 12.2 ps, z 0 b>2.5mm) of the medium (not shown in 
Fig. 6a) yield larger d because of 'edge artefacts' attributable to 
the sharp index discontinuity. When edge frames are omitted, we 
observe d<0.5 for 0.4<£ pu <0.6 uJ (compare first and second 
rows of Table 1). Second, in the last two frames for £ pu = 0.7 jij, 
we observe d>0.5. This is in part because the 'real' pulse develops 
multiple filaments (Fig. 6c, bottom) that cylindrically symmetric 
NLSE simulations fail to capture, even though FDT 



reconstructions capture a slice of them. Thus, here a shortcoming 
of the NLSE simulation, not the FDT reconstruction, contributes 
to larger d. The small index 'hole' that appears at £ Q b = 12.2ps 
(z ob = 2.5mm) also contributes, because it challenges FDT 
resolution in all three variables (z^ c \ x oh , tpu) and deviates in 
shape from the narrower (in x ob ), longer (in ) negative-index 
plasma filament that appears in the corresponding NLSE 
simulation frame. Excluding these final two frames, d < 0.5 even 
for E pu = 0.7 uJ (third row of Table 1); including them, d exceeds 
that of the phantom five-probe configuration III (see Fig. 3d) by 
only ~25%. Despite locally enhanced d, the reconstruction 
captures the essential physical feature (a transient plasma 
filament), albeit with distorted shape. Thus, direct comparison 
of Figs 6a and 7a and RMS error analysis show that the 
reconstructions are not significantly corrupted by artefacts. 

Figures 7b,c, 8 and 9 highlight several specific areas of 
quantitative agreement. First, regarding transverse profiles, Fig. 8 
compares (a) Ax OD (t OD ) (full width at half maximum) and 
(b) An(0,0,f ob ) from reconstructions (data points) with NLSE- 
simulated values (solid curves). For £ pu = 0.4, 0.5 uj> they agree 
well for all t OD ; for E pu = 0.6, 0.7 uJ, they agree up to the onset of 
multifilamentation, beyond which solid calculated curves are 
drawn dotted. Figure 7c directly compares the complete simulated 
(solid curves) and reconstructed (dashed) transverse index 
profiles near the exit plane (for E pu = 0.4, 0.5, 0.6 uJ) or at the 
onset of multifilamentation (0.7 ji J), further illustrating the close 
agreement. Second, both simulated (Fig. 7a) and reconstructed 
(Fig. 6a) longitudinal profiles evolve negligibly (Fig. 9a) and agree 
with independent autocorrelation measurements (Fig. 9b) of 
incident and transmitted pulse duration, as expected for L « L dis . 
Although longitudinal breakup of the pulse temporal envelope 
('pulse splitting') triggered by material dispersion 33 ' 40 or 
ionization 30,41 has been observed in previous fs pulse 
filamentation experiments and simulations, and is an important 
mechanism in arresting self-focusing collapse 5 , our NLSE 
simulations predict that pulse splitting would occur only for 
z ob >0.3cm (for £ pu = 0.7uJ), beyond our sample length. 
Previous observations of fs pulse splitting in glass used higher 
power pulses 30 or longer propagation lengths 33 than our 
experiment. Third, simulated (Fig. 7b) and measured (Fig. 6b) 
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Figure 7 | NLSE-simulated movies of the evolving index profile, (a) Selected 2D snapshots An(z^ cj ,x ob ) of the pump nonlinear index profile from 
NLSE simulations for the same four pump energies and same five propagation times as shown in Fig. 6a. The colour bar indicates the dimensionless 
refractive index change, (b) Spectra of the transmitted pump pulses from NLSE simulations, for direct comparison with measured spectra shown in Fig. 6b. 
(c) Solid curves: NLSE-simulated pump transverse spatial profiles at =0 near the exit plane (t ob = 12.2ps, z ob = 2.5mm) for pump energy from 
0.4 to 0.6 jiJ, and at t ob = 7.4ps (z ob = 1.5mm) for 0.7 (iJ. Dashed curves: corresponding lineouts of FDT-reconstructed pump spatial profile from 
Fig. 6a for direct comparison. 



Table 1 | RMS deviations of FDT reconstructions from NLSE 
simulations. 

t ob range z ob range E pu = 0.4 nJ 0.5 nJ 0.6 nJ 0.7 nJ 

0-14.7ps 0-3 mm 0.63 0.65 0.66 0.99 

2.5- 12.2 ps 0.5 -2.5 mm 0.37 0.49 0.50 0.95 

2.5- 7.4 ps 0.5 -1.5 mm 0.37 0.49 0.50 0.46 

FDT, frequency-domain tomography; NLSE, nonlinear Schrodinger equation; RMS, root mean 
squared. 



output pump spectra show similar trends as E pu increases. 
Specifically, the simulations reproduce the development of an 
asymmetric double-peaked spectrum for £ pu >0.6uJ, caused by 
self-phase modulation. The asymmetry in the peaks is attributable 
to slight self- steepening, which otherwise plays a minor role, as 
our sample length L is much less than the characteristic self- 
steepening length L ss = CT p Jn 2 I pu ~ 15 cm (ref. 42) for our 
conditions. 

Beyond the point at which 0.6, 0.7 uj pumps lost cylindrical 
symmetry, the calculations, although no longer in quantitative 
agreement with the data, nevertheless mirror qualitative trends in 
the FDT images. For example, for £ pu = 0.7uJ and £ ob >9.8ps 
(z ob >2.0mm), the calculated An(0, 0, f ob ) increases sharply to 
~ 0.00065 before dropping equally sharply (Fig. 8b, black dotted 
curve) due to formation of an electron-hole plasma, which 
contributes a compensating negative index change. The measured 
An (0, 0, t ob ) increases to a lower maximum, then also drops 
(Fig. 8b, black data points). Here the calculated An evolves more 
quickly than the interframe resolution; thus (like the 'dot' in the 
phantom simulations), the reconstruction 'rounds off the 
temporal or z ob dynamics, partly explaining the discrepancy. 
Absence of multifilamention in the calculation also contributes to 
the discrepancy. The FDT images themselves provide additional 
insight. In this region, the steep -walled index 'hole' appeared 
(see Fig. 6a, bottom right frame), inducing de-focusing and 
presumably pulse splitting in a longer medium 30 ' ' 40 ' 41 . 



This hole had ~ 10 Jim diameter, significantly narrower than 
the n 2 I pu envelope (see Fig. 8d, solid curve). Such a narrow hole is 
consistent with the onset of multiphoton ionization, the 
mechanism that accounts for the drop in An(0, 0, f ob ) in our 
NLSE simulation. Six photons are needed to create an electron- 
hole pair in fused silica; hence, the plasma generation rate scales 
as I 6 . Thus, plasma concentrates in the centre of the self-focused 
profile (Fig. 8d, dashed curve) and appears suddenly as 
J increases, consistent with the FDT images. Possible additional 
contributions to this feature are discussed in Methods. 

In summary, we have demonstrated a single-shot tomographic 
method for visualizing spatiotemporal dynamics of light-velocity 
objects. Compact methods to generate and detect multiple 
synchronized probe pulses were demonstrated, Radon transfor- 
mation methods were adapted to recovery of motion pictures and 
image-quality limits were elucidated through phantom simula- 
tions and laboratory experiments. FDT can image a wide range of 
nonlinear propagation phenomena, including filament formation 
in gases and the evolution of plasma wakefields. 

Methods 

Derivation of projection angle. To express the generalized Radon transform 
in the form of equation (2), it is necessary to relate object to probe coordinates. 
To do so, we first define stationary lab frame coordinates (t iy x( ,y( ,z()> where 
i = (ob,pr) and the z[ axis is parallel to the propagation direction of the object or 
probe pulse. For pulses propagating in the x-z plane, these coordinates are related 
by: z 0 b = Zp/cosfl + Xp/sinfl and x ob ' = — z pr 'smQ + x pr 'cos0. We then apply 
canonical transformations {z\ ,Xi,y b ti) = (z/ — Vit,x(,y(, f 2 ) between lab and 
pump/probe co-moving frames, where v t is the velocity of the object or probe pulse, 
to obtain 

z^ c) =z ( pr c) cos9 + x pr sin9 - (v oh - v pr cos#) t oh 
Xob= — zp? C ^ sin0 + x pr cos6 — v pr £ o bsin0 (5) 

and the inverse relations 

4r° C) = z ib OC)cos ^ _ *obSin# - (v pr - v oh cos6)t oh 
Xpr = z^sinO + x o bcos0 + v o bf o bsin0 (6) 

where we use t oh = t pr . Equation (5) immediately yield the arguments of the 
integrand in equation (2). 
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Figure 8 | Comparison of FDTdata with NLSE simulations, (a) Reconstructed radius Ax ob and (b) index shift An(z^ = 0,x ob = 0) at the centre of 

Kerr profiles An(z^ =0, x ob = 0) from Fig. 6a versus propagation time f ob . Data points are for various /f pu : red squares, 0.4 uJ; green circles, 0.5 uJ; 
blue diamonds, 0.6 uJ; black triangles, 0.7 uJ. For /f pu = 0.7uJ and f ob >8.8ps, Ax ob is not well defined because of the complex Kerr profile shape, and 
thus is not plotted. Solid curves are from NLSE calculations taking only n 2 / pu into account; the dotted portion of these curves indicates the region after 
multiple filaments formed. Dashed curves include fourth-order Kerr effect with n 4 = -3.3 x 10 ~ 28 cm 4 W~ 2 (ref. 36). Error bars in a are equal to the 
vertical axis pixel size of each movie frame in Fig. 6a and denote error in measuring x ob . Error bars in b are the maximum difference of index shift within 
interframe temporal resolution Af ob . (c) Results of z-scan measurements of Kerr lensing by the fused silica plate. Data points show far-field beam radius 
downstream of the silica plate placed at indicated time f Q b from the vacuum waist of 100 fs pulses focused to w 0 = 120 |im, / 0 = 0.6 TWcm ~ 2 in vacuum. 
Solid blue curve is calculated radius based on lowest-order Kerr lensing, dashed curve includes n 4 as above, (d) Lineout (solid blue curve) of An versus x ob 
at = 0, f ob = 12.7 ps (z ob = 2.6 mm), showing index hole near centre of the Kerr profile, compared with calculated radial distribution of plasma density 
(dashed curve) at t ob = 12.7ps (z ob = 2.6mm) due to six-photon ionization. Calculated peak plasma density overestimates actual peak density by as much 
as 2 6 because calculated / pu overestimates actual / pu at z ob = 12.7 ps (z ob = 2.6 mm) by ~ 2 (compare black dotted curve and data points in b). Dotted blue 
curve shows suppression of the Kerr profile from fourth-order Kerr effect using n 4 as above. 




Figure 9 | Direct comparison of measured and calculated longitudinal pulse profiles, (a) FDT reconstruction (blue solid line) and NLSE simulation 
(red dashed line) of the longitudinal profile, (b) Autocorrelation trace calculated from FDT-reconstructed longitudinal profile (blue solid line) and 
direct measurement (red dash line). Except for the direct autocorrelation measurement, the curves are lined out fromFigs 6 and 7 with f ob = 7.4ps, 
x ob = 0 and E pu =0.6uJ. 

To derive equation (3) relating (j) to 6, define the starting point of the streak into free space, where it is detected, these coordinates become (cz^ /v oh ,x oh ), 
axis as the common origin of Zp r oc ^ — x pr and z^ — x oh coordinate systems yielding 

at t ob = 0. The former origin moves to (z^ ,x oh ) in the latter system at t oh > 0, ^ loc ^ 
thereby defining the streak axis. Once the phase- modulated probe propagates 

0=arg[cz; b VVob + MCob], (7) 
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Similarly, the starting and end points of the streak axis have coordinates (0, 6) and 
x pr coordinates. The components of the streak 



(Zpr C \x pT ), respectively, in Zp r oc ' 



(loc) 



length vector are thus AzJ 
substituting equation (6) with z\ 



(loc) 



and Ax pr = x pr , which yields equation (8) on 
= 0 and x oh = 0 and putting L = v ob f ob . 



FDT interframe resolution limits. The length of a phase streak at probe angle 0 is 



Azjr ) + (Ax pr ) 2 , where 



Az 



(loc) 



-L cosO- 



Ax px = LsinO 



(8) 



are the two components of the streak length vector for an object propagating 
through a medium of length L. Thus, transverse profile evolution, measured using 
small 0 probes, is best resolved when Azp r oc ^ 3> Az^ c \ As Azp r oc ^ oc |v pr — v ob | at 
small 0, using pump and probe pulses with a large velocity mismatch optimizes 
transverse temporal resolution. Longitudinal profile evolution, measured using 
large 0 probes, is best resolved when Ax pr ^> Ax ob . As Axp r ocsin0, including 
probes up to 0 max = 7r/2) optimizes longitudinal zj^ resolution. As a 
numerical example relevant to the experiments, for which v ob = 0.68c, v pr = 0.66c, 



= 20 um, Ax oh = 50 um and I 



= 70°) = 14°, the temporal resolu- 



tion is <$fj b - (v pr /v ob )Az o 1 b oc) /|v ob - v pr | « 3 ps and 5% b ~ Ax ob /(v ob sin~0 ma J « 1 ps 
for transverse (T) and longitudinal (L) profile evolution, respectively. Interframe 
resolution is thus determined by probes at the extreme angles 0 and 0 max , but does 
not depend on the number or angular distribution of probes between these angles. 

Angular distribution of probes. The distribution of probe angles is a key FDT 
design feature, in addition to those discussed in the main text. For designs with 
equivalent probe number and |v pr — v ob |, an even angular distribution of probes 
with maximized (f) range yields the smallest RMS error, analogous to the limit- 
angle problem of conventional computerized tomography 25 . Absence of \cf)\ >70° 
probes in configuration III (see 'Phantom simulations' section of main text), 
illustrates this point. The absence introduces streak artefacts into the reconstructed 
profiles (Fig. 2d). The artefacts persist, and RMS error remains above that of 
configurations I, II even when probe number within \(j)\ <70° is the same as in 
those configurations (not shown). 

Experimental limits on the maximum probe number. In general, the maximum 
number of probe N is limited by the peak separation in Fig. 4c, which is related to 
the angular separation AO between probes. In the experiment, AO was limited to the 
probe divergence angle A0 min = X pr IAx oh , or ~ lOmrad for our system. On the 
other hand, the spectrometer CCD pixel size x pbL ~ 40 urn limited the maximum 
probe-reference angular separation to 0 max = M/l pr /%p ix ~2OOmrad, where M~20 
is the magnification of lens L 2 . Thus, ideally our system supports N max = 0 max / 
A# min ~20. In practice, 0 max is limited to ~ lOOmrad (and iV max to ~ 10) by the 
physical aperture and spherical aberration of lenses L x and L 2 . We chose AO 
conservatively at ~ 2A0 min to thoroughly eliminate cross-talk, yielding five probes. 

ART movie reconstruction. On the basis of equation (2), each pixel value fa of a 
2D phase streak was obtained by applying a linear projection operator P on the 3D 
tomographic movie voxels Anf. 

~~ (9) 



Following ART, we applied an iterative procedure to solve equation (9) by setting 



initial solution Am — 0 and updating this solution to fit the measured phase 
streak ij/j in the kth iteration via 



An? } =An 



(k-i) 



EmPimAn 
EmPfr, 



(k-l) 



(10) 



where X = 0.1 is the relaxation parameter for ART iterations. When the 2D 
projections ij/^\z p ? c \ x pr , y 0 ) contained In jumps, a 2D phase unwrapping pre- 
processing procedure was applied to the measured phase through the weighted 
minimum norm algorithm 43 . 

Simulations. Pump pulse propagation was simulated using the model in ref. 31, in 
which an NLSE and a coupled electron-hole plasma generation-recombination 
equation were solved numerically. 



1 dA_ i r 
v ob dt 2k 0 n 1 



iP 2 d 2 A 



+ N N 



|=^(|A| 2 ) + ^p-^ 

OL, E g T r 
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where t, = z^ loc Vv. In the NLSE, the first and second terms are responsible for trans- 
verse diffraction and longitudinal dispersion, respectively, and the third nonlinear 



term N N l = ik 0 n 2 T[(l -f R )\A\ 2 +f R R(Q (x) \A(Q\ 2 ]A - f(l + ioj 0 t c ) T~ \pA) 
— 1VpJ i^|U^A denotes instantaneous Kerr effect, delayed Raman-induced Kerr 
effect, plasma absorption/defocusing and photo -ionization loss of the laser wave. 
We also performed some simulations that included a fourth-order Kerr term n 4 I^ u 
term, as discussed in the next section. The operator T = 1 + (i/co 0 )d/dC, leads to self- 
steepening and space-time coupling. The plasma generation-recombination 
equation has taken photo -ionization, avalanche ionization and recombination 
into account, whereas the photo -ionization term Wp/(|A| 2 ) is based on 
Keldysh's formulation generalizing both multiphoton and tunnelling ionization 
mechanisms 39 . To compare simulation results with experimental data, transverse 
(longitudinal) imaging resolution Ax res ~ 5 um (Az£s /v ob ~ 70 fs) were taken into 
account by convolving the simulated index profile with a Gaussian function with 
transverse (longitudinal) dimension Ax res (Az£° C Vv ob ). Simulation parameters were 
either from our direct measurements (for example, n 2 = 2 x 10 ~ 16 W cm - 2 
from z-scan measurements) or most recent references for each parameter: 
n 4 = — 3.3 x 10 ~ 28 W cm - 2 (ref. 36); n 0 = 1.45, group velocity dispersion 
j3 2 = 36.1 fs 2 mm ~~ x , silica band gap E g = 9 eV and constants T d = 32 fs, t s = 12 fs, 
/ R = 0.18 for the Raman response function R(Q = t s (t~ 2 + t 7 2 )exp( — ^/T d ) 
sin(^/T s ) 31 ; plasma collision time t c = 1.7 fs and plasma recombination time 
T r = 170fs (ref. 44). 



High-order Kerr effect. The relative importance of plasma defocusing 12 ' 34 ' 35 and 
negative electronic Kerr effect 36-38 in arresting the collapse of self-focusing pulses 
has been debated extensively in connection with nonlinear pulse propagation in 
gaseous media. Although there has been comparatively little discussion of this issue 
for condensed media, a measured value n 4 = — 3.3 x 10 -28 cm 4 W~ 2 of the 
fourth-order Kerr coefficient of fused silica has been reported 36 . In view of this 
measurement, and the importance of nonlinear pulse propagation in fused silica in 
optical fibre communication and micro -machining, we considered the possibility 
that fourth-order Kerr effect might contribute to, or even dominate, the arrest of 
self-focusing under our conditions. There were three findings. First, our NLSE 
calculations that included fourth-order Kerr effect using the previously reported 
n 4 deviated significantly from all FDT data, even before the onset of 
multifilamentation, as shown by the dashed curves in Fig. 8a,b of the main text. 
This n 4 actually induced self-defocusing under our conditions, contrary to 
observation. RMS deviation d between the NLSE calculations and the 
reconstructions (Fig. 6a) also grew to >0.8, even in regions free of edge artefacts 
and for £ pu <0.6 uj, much larger than values in Table 1. Second, as an independent 
check of \n 4 ln 2 \> we conducted standard z-scan measurements of our silica plate's 
nonlinear index 45 . Results, exemplified by Fig. 8c of the main text (circles), fit well 
to a nonlinear lens determined solely by n 2 I pu (solid curve). Our FDT and z-scan 
results both place an upper limit of 10 ~ 29 cm 4 W ~ 2 on \n 4 \. Third, we found that a 
negative 4th-order Kerr effect, even of limited value, broadly suppressed the 
positive n 2 I pu profile, rather than producing the observed narrow index hole that 
appeared at t oh = 12.2 ps (z ob = 2.5mm) for £ pu = 0.7uJ (see Fig. 6a, fourth row, 
last panel). This is because it scales as I pu instead of Jp U , as plasma generation does. 
Thus, our results do not support an important role of a fourth-order Kerr effect. 
However, they do not rule out the possibility that a negative electronic Kerr effect 
of higher than fourth order contributes to forming this index hole. 
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